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Abstract 

This paper presents a continuous and discrete Lagrangian theory for stochastic 
Hamiltonian systems on manifolds. The main result is to derive stochastic governing 
equations for such systems from a critical point of a stochastic action. Using this result 
the paper derives Langevin-type equations for constrained mechanical systems and 
implements a stochastic analog of Lagrangian reduction. These are easy consequences 
of the fact that the stochastic action is intrinsically defined. Stochastic variational 
integrators (SVIs) are developed using a discretized stochastic variational principle. 
The paper shows that the discrete flow of an SVI is a.s. symplectic and in the presence of 
symmetry a.s. momentum-map preserving. A first-order mean-square convergent SVI 
for mechanical systems on Lie groups is introduced. As an application of the theory, 
SVIs are exhibited for multiple, randomly forced and torqued rigid-bodies interacting 
via a potential. 

1 Introduction 

Since the foundational work of Bismut [1981], the field of stochastic geometric mechanics 
is emerging in response to the demand for tools to analyze the structure of continuous and 
discrete mechanical systems with uncertainty [2; 20; 21; 28; 29; 18; 17; 23; 8]. Our specific 
goal within this context is to develop efficient, structure-preserving integrators for long- 
time simulations of stochastic Hamiltonian systems on manifolds. Our strategy is to extend 
variational integrators to this class of systems. 

Variational Integrators Variational integration theory derives integrators for mechani- 
cal systems from discrete variational principles [33; 22; 34; 27]. The theory includes discrete 
analogs of the Lagrangian, Noether's theorem, the Euler-Lagrange equations, and the Leg- 
endre transform. Variational integrators can readily incorporate holonomic constraints (via 
Lagrange multipliers) and non-conservative effects (via their virtual work) [34; 27]. Alto- 
gether, this description of mechanics stands as a self-contained theory of mechanics akin to 
Hamiltonian, Lagrangian or Newtonian mechanics. 

One of the distinguishing features of variational integrators is their ability to compute 
statistical properties of mechanical systems, such as in computing Poincare sections, the 
instantaneous temperature of a system, etc. For example, as a consequence of their varia- 
tional construction, variational integrators are symplectic [9; 32; 12]. A single-step integrator 
applied to a mechanical system is called symplectic if the discrete flow map it defines ex- 
actly preserves the continuous symplectic 2-form and is otherwise called standard. Using 

*Applied & Computational Mathematics (ACM), Caltech, Pasadena, CA 91125 
(nawaf @acm . caltech . edu) . 

t Applied & Computational Mathematics (ACM) and Control &; Dynamical Systems (CDS), Caltech, 
Pasadena, CA 91125 (owhadiOacm.caltech.edu). 



1 Introduction 2 



backward error analysis one can show that symplectic integrators applied to Hamiltonian 
systems nearly preserve the energy of the continuous mechanical system for exponentially 
long periods of time and that the modified equations are also Hamiltonian (for detailed 
exposition see [13]). Standard integrators often introduce spurious dynamics in long-time 
simulations, e.g., artificially corrupt chaotic invariant sets as illustrated in a computation 
of a Poincare section of an underwater vehicle in Fig. 1.1 using a fourth-order accurate 
Runge-Kutta (RK4) method and a variational Euler ( VE) method [4] . 




(a) RK4 (b) RK4 (c) VE (d) VE 

h=0.025 h=0.05 h=0.025 h=0.05 



Figure 1.1: Underwater Vehicle Dynamics [4]. This figure shows a computation of Poincare 
sections using a second-order accurate variational Euler integrator (VE) as compared to fourth- 
order accurate Runge-Kutta (RK4). Both methods agree with the benchmark at the finer stepsize 
h — 0.025. However, at the coarser stepsize h = 0.05, RK4 corrupts chaotic invariant sets while the 
lower-order accurate VE method preserves the structure of the benchmark. 

In addition to correctly computing chaotic invariant sets and long-time excellent en- 
ergy behavior, evidence is mounting that variational integrators correctly compute other 
statistical quantities in long-time simulations. For example, in a simulation of a coupled 
spring-mass lattice, Lew et al. found that variational integrators correctly compute the time- 
averaged instantaneous temperature (mean kinetic energy over all particles) over long-time 
intervals, whereas standard methods (even a higher-order accurate one) exhibit a artificial 
drift in this statistical quantity. See Fig. 1.2 and [19]. These structure-preserving properties 
of variational integrators motivate their extension to stochastic Hamiltonian systems. 
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Figure 1.2: Coupled Spring— Mass Lattice [19]. This figure shows several plots of the time- 
averaged instantaneous temperature as a function of time using a first-order accurate variational 
integrator (VII) and RK4. The mechanical system is a coupled spring-mass lattice. Observe the 
drift in the computed temperature by RK4 as compared to VII. 



Main Results In his foundational work, Bismut showed that the stochastic flow of cer- 
tain randomly perturbed Hamiltonian systems with configuration space K. n extremize a 
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stochastic action. He called such systems stochastic Hamiltonian systems, and used 
this property to prove symplecticity and extend Noether's theorem to such systems [2]. 
Mean-square symplectic integrators for stochastic Hamiltonian systems with configuration 
space W 1 and driven by Wiener processes have been developed [28; 29]. 

Bismut's work was further enriched and generalized to manifolds by recent work [18; 
17]. Lazaro-Cami and Ortega showed that stochastic Hamiltonian systems on manifolds 
extremize a stochastic action defined on the space of manifold- valued scmimartingales [18]. 
Moreover, they performed reduction of stochastic Hamiltonian systems on the cotangent 
bundle of a Lie group to obtain stochastic Lie-Poisson equations [17]. However, as far as we 
can tell the converse to Bismut's original theorem, namely a critical point of a stochastic 
action satisfies stochastic Hamilton's equations, has not been proven. In fact, as pointed out 
by Lazaro-Cami and Ortega, a counterexample can be constructed to prove the converse of 
Bismut's theorem is not true for a certain choice of stochastic action [18]. 

In this paper we restrict our attention to stochastic Hamiltonian systems driven by 
Wiener processes and assume that the space of admissible curves in configuration space are 
of class C . From the viewpoint of randomly perturbed mechanical systems, this latter 
restriction is reasonable since random effects often do not appear in the kinematic equation, 
but rather the balance of momentum equation as white noise forces and torques. It should 
be mentioned that the ideas in this paper can be readily extended to stochastic Hamiltonian 
systems driven by more general semimartingalcs, but for the sake of clarity we restrict to 
Wiener processes. Within this context the results of the paper are as follows: 

• For a class of mechanical systems whose configuration space is a paracompact manifold 
and which is subjected to multiplicative white noise forces and torques, the paper 
proves a.s. a curve satisfies stochastic Hamilton's equations if and only if it extremizes 
a stochastic action. This theorem is the main result of the paper. 

• The paper derives governing SDEs for stochastic Hamiltonian systems with holo- 
nomic constraints using a constrained stochastic variational principle and for stochas- 
tic Hamiltonian systems with nonconservative effects in the drift using a Lagrange- 
d'Alembert-type principle (for deterministic treatments see [24]). The paper performs 
Lagrangian reduction for stochastic Hamiltonian systems whose configuration space 
is a Lie group and provides stochastic Euler-Poincare/Lie-Poisson equations for such 
systems (for deterministic treatment see [24] ) . These are easy consequences of the fact 
that the stochastic action is intrinsically defined. 

• The paper shows how to discretize variational principles to obtain stochastic varia- 
tional integrators (SVIs), stochastic RATTLE-type integrators for constrained stochas- 
tic Hamiltonian systems, and stochastic Euler-Poincare/Lie-Poisson integrators for 
stochastic Hamiltonian systems on Lie Groups (for deterministic treatments see [31; 
34; 25; 13]). In addition, the paper describes how to derive quasi-symplectic methods 
for rigid-body-type systems at uniform temperature. 

Organization of the Paper Sufficient conditions for existence, uniqueness, and mean- 
squared differentiability of stochastic flows on manifolds are recalled in §2. In §3 we extend 
the Hamilton-Pontryagin (HP) principle to the stochastic setting to prove a class of mechan- 
ical systems with multiplicative noise appearing as forces and torques possess a variational 
structure. It should be emphasized (and it is explained in the section) the mechanical 
system could evolve on a nonlinear configuration space and involve holonomic constraints 
or nonconservative effects in the drift. The HP viewpoint is adopted since it unifies the 
Hamiltonian and Lagrangian descriptions of the system. By left-trivializing this principle, 
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we also show how to perform Lagrangian reduction in this stochastic setting for stochastic 
rigid-body-type systems. 

In §4, SVIs are derived from an abstract discrete Lagrangian and the structure of the 
resulting discrete flow map is analyzed. In §5, we concretely show how to design a single-step, 
stochastic variational Euler integrator for mechanical systems whose configuration space is 
a Lie group using a simple stochastic discrete HP principle. If the configuration space is 
the translation group, the governing SDEs and integrators are in one-to-one correspondence 
with so-called stochastic Hamiltonian systems and mean-square symplectic integrators [2; 
28; 29]; and with the addition of dissipation in the drift term, in one-to-one correspondence 
with Langevin equations and so-called quasi-symplectic methods for such systems. These 
mean-square symplectic integrators have been numerically tested and shown to possess 
excellent properties for long-time simulations of mechanical systems governed by Langevin- 
type equations; see [30]. 

Our own simulations confirm those findings, and will be discussed in [5; 6; 7]. A sample 
of such results is provided in Fig. 1.3. It compares an SVI to standard, presumably non- 
variational methods on a ballistic pendulum at uniform temperature [5]. The figure shows 
that an SVI correctly computes the temperature of the system (defined as the mean of 
the instantaneous temperature with respect to realizations) whereas explicit and implicit 
Euler-Maruyama schemes (EEM, IEM) do not. All of these methods are first order strongly 
convergent and driven by the same realization of Wiener processes. This computation 
suggests that an SVI has favorable long-time energy properties whereas EEM and IEM 
artificially heat and cool the system respectively. On the other hand, this paper focuses on 
SVI theory and the structure-preserving properties of SVIs. 

In §6, as an application we explain how one can add multiplicative white noise forces 
and torques to multiple rigid bodies in a fashion that preserves variational structure. With 
the addition of dissipation these become Langevin-type equations. An SVI is provided for 
such systems. Using an Ito- Taylor expansion, it is easy to check that the method is first 
order in the mean-squared sense [15]. 




(a)/i = 0.1 (b) h = 0.05 



Figure 1.3: Ballistic Pendulum at Uniform Temperature [5]. Plots of the mean instan- 
taneous temperature (kinetic energy) of a ballistic pendulum computed using an SVI, explicit 
Euler-Maruyama (EEM), and implicit Euler Maruyama (IEM) for time-steps h as indicated. The 
correct temperature is indicated by the solid line. Observe that the EEM and IEM schemes arti- 
ficially heat and cool the system respectively. A key feature of the ballistic pendulum is that the 
diffusion and drift matrices associated to the momentums are degenerate, yet the system is still at 
uniform temperature. 
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2 Stochastic Flows on Manifolds 

Some standard results on flows of SDEs on manifolds are reviewed here for the reader's 
convenience. The reader is referred to the following textbooks on the subject for more 
detailed exposition [10; 11; 14; 16]. This section parallels the treatment of deterministic 
flows on manifolds found in Chapter 4 of [1]. 

We start by introducing notation for deterministic vector fields on manifolds which are 
an important component of SDEs on manifolds. Let M be an n-manifold. Recall, a vector 
field on M is a section of the tangent bundle TM of M. The set of all C k vector fields on 
M is denoted by X k {M). 

The notion of a probability space is introduced in order to extend the definition of a 
dynamical system to incorporate noise. A stochastic dynamical system consists of a base 
flow on the probability space which propagates the noise and a stochastic flow on M which 
depends on the noise. 

Definition 2.1 (Stochastic Dynamical System). A stochastic dynamical system consists 
of a base flow on a probability space (f2, J 7 , P) and a stochastic flow on a manifold M. The 
base flow is a ¥ -preserving map, 6 : M. x Q — ► Q, which satisfies: 

1. 9q = idn : £1 — > Q, is the identity on fl 

2. for all s, t £ K the group property: 9 S o Q t = 8 s+t 

Given times < r < s < t, the stochastic flow on M is a map (pt s : £1 X M — ► M such 
that 

1. for a. a. lo G Q, the map (s,t,uj,x) i— > ip ttS (u;)x is continuous in s, t and x 

2. ip SjS (u>) = idu '■ M — > M , is the identity map on M for all s £ K 

3. ip satisfies the cocycle property: 



This paper is concerned with stochastic dynamical systems that come from stochas- 
tic laws of motion; that is, ones whose stochastic flows define solutions of SDEs. The 
Stratonovich definition of stochastic integrals is adopted to extend SDEs to a manifold M. 
The main advantage of the Stratonovich approach being that the chain rule holds for the 
Stratonovich differential. Consider a manifold M modelled on a Banach space E, and vector 
fields Xi € X k (M) for i = 0, m. Let Tt be a nondecreasing family of er-subalgebras of T , 
and (Wj(t, cj), J-'t), i — 1, ...,m, be independent Wiener processes for < t < T,. In terms 
of these objects, the Stratonovich SDE the paper considers takes the form, 



Xq is referred to as the drift vector field and Xi, i — l,...,m, are the diffusion vector 
fields. A Stratonovich integral curve of (2.1) is a C°-map c(-,ui) : [0,T] — > M which 
satisfies: 




(2.1) 




Uniqueness of solutions to (2.1) will be defined in the pathwise sense. 
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Definition 2.2 (Pathwise Uniqueness). Let c be a Stratonovich integral curve of (2.1). 
Pathwise uniqueness of c means that ifc:I^Mis also a solution to (2.1) on the same 
filtered probability space with the same Brownian motion and initial random variable, then 

P(c(t,uj) = c(t,u), V <t <T) = 1. 

Differentiability of the flow map on M will be defined in the mean-squared sense. In the 
following we define mean-squared derivatives on the model space E with the understanding 
that this notion can be extended to M using a local representative of the flow map. 

Definition 2.3 (Mean-Squared Derivative). The mean squared norm of f : E x f2 — > E is 

given by: 

H/(.x,^)|| = (E(|/(x,^)| 2 )) 1/2 . 

Using this norm one can define the derivative of f in the standard way, i.e., f is mean 
squared differentiable at a £ E if there is a bounded linear map Df(a,ui) : E — > E that 
satisfies, 

,. \\f(a + 6,Lj)-f(a,Lj)-Df(a,u;)-6\\ 

J lm iTxTi * °- 

5^0 \\S\\ 



As is standard the explicit dependence on the point to £ ft will usually be suppressed. 
With these definitions one can state the following key, but standard theorem [10; 11; 14; 16]. 

Theorem 2.4 (Existence, Uniqueness, and Smoothness). Let M be a manifold with 
model space E. Suppose £ X k (M), i = 0, ...,m and k > 1, are uniformly Lipschitz, and 
measureable with respect to x £ M. Let L — [0,T]. Then, 

1. For each u £ M , there is a.s. a C Q -curve c : I — » M such that c(0) = u and c satisfies 
(2.1) for all t £ I. This curve c : I — > M is called a maximal solution. 

2. c is pathwise unique. 

3. There is a.s. a mapping F : I X M — > M such that the curve Cu : I — * M defined 
by c u (t) = F t (u) is a curve satisfying (2.1) for all t € I. Moreover, F is C k in the 
mean-squared sense in u and C° in t. 



3 Stochastic HP Mechanics 

In this section a variational principle is introduced for a class of stochastic Hamiltonian 
systems on manifolds. The action is simply the classical action stochastically perturbed 
by stochastic integrals. The key feature of this principle is that one can recover stochastic 
Hamilton's equations for these systems. Roughly speaking, this is accomplished by means of 
taking variations of this action within the space of curves only (not the probability space), 
and imposing this "partial differential" of the action to be zero. 

Stochastic HP Principle Let (f2, P) be a probability space. Fix an interval [a, b] C R 
and let T t be a nondecreasing family of er-subalgebras of T, and (Wi(t), !Ft), i = 1, m, be 
independent Wiener processes for t £ [a, b]. Consider a mechanical system with Lagrangian 
C : TQ — > K and a paracompact, configuration manifold Q. 

The paper adopts an HP viewpoint to describe this mechanical system with random 
perturbations. The HP principle unifies the Hamiltonian and Lagrangian descriptions of a 
mechanical system [35; 36; 3; 4]. The classical HP action integral will be perturbed using 
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stochastic potentials: 7$ : Q — ► R for i = l,...,m. Roughly speaking, in the stochastic 
context the HP principle states the following critical point condition on TQ © T*Q, 



C{q,v)dt + ^ ll {q)odW t + 
=1 



"■ ft ~ 



= 0, 



where (q(t),v(t),p(t)) G TQ © T*Q are varied arbitrarily and independently with endpoint 
conditions q(a) and q(b) fixed. This principle builds in a Legendre transform, stochastic 
Hamilton's equations and stochastic Euler-Lagrange equations. The action integral in the 
above principle, consists of two Lebesgue integrals with respect to the Lebesgue measure 
dt and m Stratonovich stochastic integrals with respect to Wiener processes. This action 
is random, i.e., for every sample point to G f2 one will obtain a different, time-dependent 
Lagrangian system. However, each system possesses a variational structure which we will 
make precise in this section. For a deterministic treatment of time-dependent continuous 
and discrete Lagrangian systems the reader is referred to [27]. 

Definition 3.1. The Pontryagin bundle is defined as the Whitney sum TQ(&T*Q. Fixing 
the interval [a, b], define path space as 

C{[a,b],q 1 ,q 2 )={(q,v,p)eC°{[a,b],TQ®T*Q) | q G C\[a, b], Q), q(a) = q 1} q(b) = q 2 }. 
Let (& : x C([a, b], q\, q 2 ) — > M. denote the stochastic HP action integral. 



<8(u,q,v,p) = / 

J a 



£{q,v)dt+J2^ i {q)odWi{t,cj) 



The Pontryagin bundle is a vector bundle over Q whose fiber at q S Q is the vector 
space T q Q © T*Q. The HP path space is a smooth infinite-dimensional manifold. One 
can show that its tangent space at c = (q, v,p) £ C([a, 6], q±, q 2 ) consists of maps w = 
(q,v,p,5q,6v,6p) G C°([a,b],T(TQ © T*Q)) such that 5q{a) = 5q(b) = and q,Sq are of 
class C 1 . Let (q,v,p)(-,e) G C(qi,q 2 , [a,b]) denote a one-parameter family of curves in C 
that is differentiable with respect to e. Define the differential of 25 as, 



d 

d(5 • (Sq,Sv,6p) := — <S(u,q(t,e),v(t,e),p(t,e))\ e=0 



where 
Sq(t) -- 



d d d 

-?(*,e)le=o. 6q(a)=6q(b) = 0, 5v(t) = — v(t, e)| e=0 , Sp(t) = — p(t, e)| e=0 . 



In terms of this differential, one can state the following critical point condition. 

Theorem 3.2 (Stochastic Variational Principle of Hamilton-Pontryagin). Let £ : 

TQ — > K be a Lagrangian on TQ of class C 2 with respect to q and v and with globally 
Lipschitz first derivatives with respect to q and v. Let 7, : Q — * K be of class C 2 and with 
globally Lipschitz first derivatives for i = 1, m. Then almost surely, 

A curve c— {q,v,p) G C([a, b], q%, q 2 ) satisfies the stochastic HP equations: 

vdt, 




^dt 

dq UL 
dC 



m d-fj 
dq 



Em 
i=l 



■dWi. 



(3.1) 



if and only if it is a critical point of the function : x C([a, b],qi, q 2 ) 
that is, d<S(w,c(w)) = 0. 
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Proof. Let us first prove almost surely: a critical point of the function satisfies (3.1). The 
differential of the HP action integral is given by, 



d®(uj,c)-(5q,Sv,Sp) 
dC 



dv 



Svds 



Sp, 



dC 

dq 
dq 



dt 



»=i * 

— v ) ds + ( p, S 



■ Sq o dW, 
dq 



dt 



Sv ) ds 



One can use a dominated convergence argument to show that differentiation and stochastic 
integration commute in the above stochastic integrals, as 7$ are of class C 2 for i = 1, ...,m, 
and the curves are continuous. Consider the term involving Sp. Since Sp is arbitrary and the 
integrand is continuous, the kinematic constraint holds: dq/dt = v. 1 Similarly, the Legendre 
transform is obtained from the Sv term: dC /dv = p. 

Collecting the variations with respect to Sq in the differential gives, 



dC 
dq 



Sqds + (p,S 



.dq 
> — 
dt 



m „ 

V — 



Sq o dWi 



The first two terms are standard Lebesgue integrals and the last m terms are Stratonovich 
stochastic integrals. Motivated by Riemann-Stieljes integration by parts the following defi- 
nition is introduced: 

Definition 3.3. Given /1 G C°([a,b],T*Q) and f 2 € C 1 ([a,b],TQ) , define 



(dfi,f 2 ) := (/i,/a)|* 



h df d l) [«■!>]■ 



Using this definition and the boundary conditions Sq(a) — 5q{b) = 0, the following 
function / : C(gi, <?2, [a, b]) x C 1 ([a, b], TQ) — > K is introduced, 



i(q,v,p, f) 



dC 
dq 



ds + £ 7T dWi - d P)-f 



so that, 



I(q,v,p,Sq) 



dC 
dq 



Sqds + (pJ^) ds + ' 6 1° dWi 



In the following it is shown that if I(q, v,p, f) — for arbitrary / of class C 1 then (q, v,p) 
satisfy (3.1). 

Let {U a ,g a } be a partition of unity on TQ © T*Q. Expand / in terms of this partition 
of unity, 



g a (q,v,p) + -a 1 dWi -dp \ ■ f 



dq 



Since the curves (g, v,p) are compactly supported, only a finite number of the g a are nonzero. 
For each g a nonzero, the terms in the integral can be expressed in local coordinates. Observe 
that since dq = vdt, the Stratonovish-Ito conversion formula implies that, 



6 97i A JTI/ 

g<x-^z ■ oq o dWi 



dq 



9a— ■ SqdWi 



dq 



1 This follows from the basic lemma that if f,g £ C°([a,b],M) and g is arbitrary then f b f(t)g(t)dt 
f(t) = V t 6 [0,6]. 
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for i = 1, m. 

We will select / to single out the jth-component of the covector field in /. Introduce the 
following function h : R — > R for this purpose: 



h(t) = 2 



t t 2 



Observe that h(0) = 0, h(e) = 1, and h'(e) = 0. Let {e^} be a basis for the model space of 
Q. Now fix j, and define / e 6 C' 1 ([a, b],TQ) in local coordinates as follows: 



fe(s) 



h(s — a)ej if a < s < a + e, 

: if a + e<s<t — e, 

— s)e_,- if t — e < s < t, 

if t < s < b. 



Introduce the following label to simplify subsequent calculations 
A(s) 



dq 

In terms of A(s), one can write 



' ^ v(s))ds + ^(<7(*)) dWi(s) - dp(s) 



h(s - a)g a (s)A(s) + / g a (s)A(s) + / fi,(t - s)g a (s)A(s) 



The proof shows in the mean squared norm (cf. definition 2.3), 



]xml(q,v,p,f e ) = V / 3 Q A(s) =: /* 



(3.2) 



Using this result and the Borel-Cantclli lemma, one can deduce there exists {e„} that 
converges to such that I(q,v,p, f 6n ) a.s. converges to /*. It follows that I* — almost 
surely. 

We proceed to prove (3.2). Since (a + b) 2 < 2a 2 + 2b 2 , 



^2 / SaA( s ) - I(q> v ,Pi fe) 
™ Jo. 

I + (1 - h(s - a))g a A(s) + f (1 - h(t - s))g a A(s) 

a Ja Jt-e 

E P + \l - Hs ~ a))g a A(s) 



< 2 



W (1 - h(t - s))g a A(s) 

„. Jt-e 



We will only show how to bound the first term since bounding the second term is very similar. 
By continuity of {q, v,p), one can pick e small enough so that the support of (q, v,p) lies in 
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E 



a+e 



a single chart. Therefore, 

(1 - h(s - a))g a A(s) 

'ac 



< 3 



(1 — h(s — a))^^ds 



i=i 



a+e 



(l-h(s-a))^j dWi 



a+e 



(1 — h(s — a))dp ■ ej 

Since is continuous on s G [a, a + e], the first term can be bounded, 

M 2 e 2 



(1 — his — a))— — rds 



< 



Similarly, by the Ito isometry and since g^j is continuous on s G [a, a + e], the second 
m terms can similarly be bounded, e.g., the ith Stratonovich integral can be bounded as 
follows, 



(l-h{s-a))^ dWi 



(•a-\-c 



E 



(1 - h(s-a)) 



ch/i 
dqi 



ds < 



M 2 t 



Using definition 3.3 and the integral mean value theorem, the final term can be bounded as 
well: 



a+e 



(1 — h(s — a))dp ■ e 3 



(1 _ h{s _ a))p 3 ( S )\lZ7 £ + I pMh'is - a)d 
-Pj(a) +pj(a + ce)\\ 2 



where < c < 1 is a real constant. Since pj is of class C°, as e — * this term vanishes. 
Since j is arbitrary we have proved (3.2). Therefore, almost surely: if c is a critical point 
of 25 then d(3(w, c) • w = for all w € T c C([a, b], <7i,(fe), and hence, c satisfies the stochastic 
HP equations. 

On the other hand, almost surely: if c satisfies (3.1), then it is a critical point of 0. 
This direction is easy to confirm, since as a solution to the stochastic HP equations c 
is a measureable diffusion process. In fact, this direction is similar to the one Bismut 
originally established, namely that the solution of stochastic Hamilton's equations extremize 
an action function; albeit the stochastic action used by Bismut has a different domain than 
the stochastic action used in this proof [2] . ■ 

Equations (3.1) are a stochastic differential algebraic system of equations. Assuming one 
can eliminate v using the Lcgcndrc transform, these equations can be viewed as a Cauchy 
problem. This paper is primarily concerned with forces or torques that appear as white 
noise in the balance of momentum equations, which explains the choice 7^ = ji(q). Observe 
that by the Ito-Stratonovich conversion formula, the Ito modification to the drift is equal 
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to zero, and hence, (3.1) can be written in Ito form as: 

dq = vdt, 

d P = ^dt+Y d ^dw i , 



P dv ' 

In what follows structure-preserving properties of the flow map defined by the maximal 
solution of these equations over [a, b] will be investigated. First, observe that because of 
the smoothness conditions assumed in theorem 3.2, a solution almost surely exists and is 
pathwise unique on [a, b] by the results in §2. When 7, is constant for i = 1, to, the reader 
is referred to the following texts for deterministic treatments of symplecticity, momentum 
map preservation, and holonomically constrained mechanical systems: [24; 27]. 

Symplecticity Assuming that one can eliminate v using the Legendre transform, the 
stochastic HP equations define a stochastic flow on the symplectic manifold (T*Q, n) where 
K is the canonical symplectic form [24]. We denote this flow by F t : T*Q — ► T*Q. With 
this assumption and in a more general context, Bismut extended the variational proof of 
symplecticity and Noether's theorem to stochastic Hamiltonian systems [2]. In fact, one 
does not need to prove both directions of (3.2) to perform this extension. These proofs are 
repeated here in the context of stochastic Hamiltonian systems driven by Wiener processes 
for the reader's convenience and for completeness. 

The variational proof of symplecticity will be used to show this flow preserves k [24]. 
By theorem (2.4), assuming the Lagrangian is sufficiently smooth then F t : T*Q — > T*Q 
is a.s. differentiablc. Recall the derivative is defined in the mean squared sense (cf. defini- 
tion 2.3). The idea of the proof is to restrict 25 to the space of pathwise unique solutions, 
i.e., define © = 25| so i ut j ons - On the same filtered probability space with the same Brownian 
motion, this solution space can be identified with the set of initial conditions, i.e., this re- 
stricted action can be expressed as, <8(w, •) : T*Q — > R. For each initial condition by (2.4), 
there exists a pathwise unique solution almost surely One then computes d©: 

d®(u,q(a),p(a)) ■ (Sq(a),5p(a)) = 

b[/ % dt + P i ^q 0dWi - d p)- 6 'i + S P-(ft- V ) ds+ {^-p)- Svdt 
+ {P,Sq)\ b a . 

The integral in the above vanishes since is restricted to solution space. The boundary 
terms define in local coordinates the canonical 1-form O on T*Q. Computing d 2 (3 gives 
conservation of the canonical symplectic form. 

Theorem 3.4 (Conservation of stochastic symplectic form). Provided that one can 
eliminate v using the Legendre transform, the flow of (3. 1 ) preserves the canonical symplectic 
form almost surely. 



Noether's Theorem In what follows we review for completeness Bismut's extension of 
Noether's theorem [2]. Let G be a Lie group with Lie algebra g. The left action of the Lie 
group on Q is denoted <!> : G x Q — > Q. The cotangent lift of this action is likewise denoted 
$ T * Q : G x T*Q -> T*Q: 

® T 'Q(h,q,p) = (*(fc >g ),(P 9 *(M))-T-P) • 
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The corresponding infinitesimal generators are £® : Q — ► T*Q and £ T ^ : T*Q — > T(T*Q) 
and by definition , 

£ Q (<Z,p) - ^ [*(exp(aO, ?)Uo - £ T * Q (?,P) = ^ [<J> T * Q (exp(^), <?,?)] s=o ■ 
This action gives rise to the following momentum map J : T* Q —> q* 

J(q,p)< = (pA Q (g,p)), 

The following conservation law follows if © is infinitesimally invariant with respect to the G- 
action. We remark in passing that infinitesimal invariance of © follows from left-invariance 
of the stochastic HP action with respect to the G-action. 

Theorem 3.5 (Stochastic Noether's Theorem). Let G be a Lie group. Lf(& is infinites- 
imally symmetric with respect to the left (or right) action of G, then the corresponding mo- 
mentum map is conserved a.s., i.e., J = (p, £q(<?)), is a conserved quantity under the flow 
of (3.1). 

Proof. Consider the differential of (S in the direction of £ T ® , 

d®(u,q(a),p(a))-e* Q (q(a),p(a)) = (p, £ Q p) > | * . 

Moreover infinitesimal symmetry implies that 

d®-e' Q (q(a),p(a)) = 0, 
=> J(q(b),p(b)) ■ £ - J(q(a),p(a)) •£ = (), 

and hence J is conserved under the flow since £ is arbitrary. ■ 

Holonomic Constraints The following results will require the proof of the converse of 
theorem 3.2. The setting in this part is an n-manifold Q and a stochastic Hamiltonian system 
with holonomic constraint. To be specific, suppose that the motion of the mechanical system 
is constrained to a submanifold S C Q defined as S = g _1 (0) where g : Q — > ]R fc , k < n, g 
is smooth, and is a regular value of g. 

As opposed to using generalized coordinates on TS, we wish to describe the mechanical 
system using constrained coordinates on TQ and introduce Lagrange multipliers to enforce 
the constraint. However, because of the stochastic component of the action, the standard 
Lagrange multiplier theorem will not apply directly and one cannot introduce Lagrange 
multipliers in the standard way. Instead, we will introduce the Lagrange multiplier using 
definition 3.3. 

In particular, consider the following constrained variational principle, 

6{<d+ [ (d\,g)) =0, 



where using definition 3.3 



(d\,g) := (A,5)|*-Jf (A, ^) dt, te[a,b}. 

In this case X(t) is a Lagrange multiplier dual to g(q(t)) for t E [a,b], and we assume is of 
class C . The corresponding equations of motion are obtained in a similar fashion as (3.1). 
In particular, one can prove the following equivalence: 



dg 
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Theorem 3.6 (Constrained, Stochastic HP Principle). Given a stochastic Hamilto- 
nian system with Lagrangian C : TQ — > R such that d 2 C/dv 2 is invertible, stochastic poten- 
tials 7j- : Q — > K for j = 1, ...,m, and holonomic constraint g : Q —> M. k with S — g _1 (0). 
Set £ = C\ TS and 71? = jj\ s for j = 1, ...,m. Fix (71,(72 G <5 afid let i : S ^ Q be the 
inclusion mapping. Then the following are equivalent: 

(i) Stochastic HP principle for C s and 7J, j = 1, m, on TS T* S (cf. thm. 3.2). 

(ii) There exists A S C°([a, swc/i £/iai z = (q,v,p) £ C([a, b], i(qi), ^(32)) A 

extremize the augmented action A) = 25(z) + (A, ^z)) where $(z)(t) = dg / dt{q(t)) 
and (\^(z))=^(\(s)^(z)( S ))d S . 

(iii) There exists A 6 C°([a, suc/i £/ia£ z = (q,v,p) £ C([a,b],i(qi),i(q2)) and A 
satisfy the constrained, stochastic HP equations 

= vdt, 

= ^(«.«), 
= 0. 




From this equivalence it follows that the flow of (3.3) is mean-square symplectic. For a 
proof of this theorem the reader is referred to [6] . 

Nonconservative Effects Nonconservative effects are incorporated by considering the 
Lagrange-d'Alembert-Pontryagin principle. In this principle the effect of a noncon- 
servative force appears as virtual work. Consider a force field F : TQ — > T*Q. Then the 
Lagrange-d'Alembert-Pontryagin principle is given by 



C{q,v)dt + J^-fi(q)odWi 



i=l 



p. 



dq 
dt 



v ) dt 



+ / F(q,v)-Sqdt = 0, 



where the variations of the first term are the usual ones (vanishing at the endpoints) . This 
principle provides a simple way to add dissipative effects into the drift which, e.g., appear 
in the standard Langevin equations for particles. 



Lagrangian Reduction For background and exposition on Lagrangian reduction in the 
deterministic setting the reader is referred to [26; 24]. Suppose that Q is a Lie group G 
with Lie algebra q. In this context one can define a left-trivialized Lagrangian by using the 
left-action of G to left-trivialize £. One does this by transforming a point (g(t), v(t)) £ TG 
to (g(t),£,(t)) £ G x g via the relation between the velocity at g(t) £ G and the body 
angular velocity at e £ G given by: £(£) = g(t)~ 1 v(t) £ g. Denote by / : G x g — > R 
the deterministic left-trivialized Lagrangian defined as l(g(t),^(t)) = C(g(t), g(t)£(t)). The 
variational principle associated with I is the left-trivialized HP principle which can be 
written as: 

-idg 
dt 



l(g,0dt- 



•5>i(fl)°<Wi 



= 



In this principle the Lagrange multiplier /i £ g* is the body angular momentum. For 
more details on the geometry of this principle in the deterministic setting the reader is 
referred to [4; 3]. The resulting equations are obtained by taking arbitrary variations with 
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fixed endpoints on g. For a function U : G 

U g G g*, as 

'dU 



u a -v 



dg 



In terms of U g , one can write the stochastic left-trivialized HP equations: 



dg 
dt 



dfi = Igdt + 2j(7i)s ° <Wi(w, £), 



i=l 



define its left-trivialized differential, 

(3.4) 

(3.5) 
(3.6) 

(3.7) 



By eliminating £ using (3.7), one obtains a SDE on G x g*. The kinematic constraint in 
this context is referred to as the reconstruction equation. We summarize this section with 
the following theorem. 

Theorem 3.7 (Stochastic Left-Trivialized Hamilton-Pontryagin Principle). Con- 
sider a mechanical system on a Lie group G with left-trivialized Lagrangian I : G x q —>■ W. 
and stochastic potentials ji : G — > K for i = 1, ...,m. Let s denote the left-trivialzied action 
given by, 

s = J" l(g,£)dt + H{g) o dW t + ^9- ld f t - dt 

Almost surely, the stochastic HP principle on a Lie group (cf. thm. 3.2), is equivalent to the 
stochastic left-trivialized HP principle: 

5s = 0, 

where the curves 

g{t) €G,£(t) €Q,u(t) eg*, te[a,b] 

can be varied arbitrarily with Sg(a) = 5g(b) = 0. A curve is a critical point of the left- 
trivialized action if and only if it satisfies the left-trivialized HP equations given by (3.5)- 
(3.7). 



4 Stochastic Variational Integrators 

The standard approach of deriving variational integrators is extended to the stochastic 
context in this section; see, e.g., [27]. The cornerstone of variational integration theory 
is the discrete Lagrangian. In this section we develop and analyze integrators from an 
abstract discrete Lagrangian which takes values on the configuration space squared. In 
the subsequent sections, discrete Lagrangians will be specified and the associated time- 
integrators analyzed from the HP viewpoint. Let [a, b] and N be given, and define the fixed 
step size h = (b— a) /N and tk = hk + a, k = 0, N. 

Definition 4.1. Consider a mechanical system with given discrete Lagrangian Cd '■ QxQ — > 
R. Let 8 t : f2 — > Q, t £ [a, b] denote the base flow on the probability space (cf. definition 2.1). 
Let the approximant to the Stratonovich integrals be denoted by: 

m rt k+ i 

B d (t k ,q k ,q k+ i,uj)hK^ / 7i ((?(*)) ° dWi(t, tk+i (u)) 

z=l 
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The associated stochastic discrete Lagrangian LrfiRxflxQxQ^Rw defined as: 

L d {t k ,u,q k ,qk+i) = £(<?fc,<?fc+i) + B d (t k ,q k ,q k +i,uj). 
Fixing the interval [a,b], define the discrete path space as 

C d (Q) = {q d - {t k }^ = - <?}■ 
Let (3 d : ft x C d (Q) — > K denote the stochastic action sum, 

N-l 



<8 d (uJ,q d ) = ^2 L d {t k ,u),q k ,q k+1 )h. 



The discrete stochastic Hamilton's principle states that the path the mechanical 
system takes in C d is one that extremizes & d (ui, •) subject to fixed endpoint conditions 
go = q(a) and qw = q(b). By discrete integration by parts (re-indexing), 



N-l 



d& d (uJ,q d ) ■ {Sq k } = ^2 (D 3 L d (t k ,uj,q k ,q k+1 ) + D A L d (t k -i, u, q k -i, q k )) ■ Sq k 
k=l 

+ D 3 L d (to,u, q ,qi) ■ Sq + D 4 L d (t N -i, w, gjv-i, Qn) -Sq N . 
Using the endpoint conditions 5qo = Sqjy = one obtains: 

N-l 

d<6 d (w, q d ) ■ {Sq k } = ^2 (D 3 L d (t k , u, q k , q k+ i) + D 4 L d (t k -i,u, qk-i,qk)) ■ Sqk- 
k=i 

Stationarity of this action sum implies the following stochastic discrete Euler- Lagrange 
equations: 

D 3 L d (t k ,uJ,q k ,qk+i) + D 4 L d (t k -i,u>,q k -i,q k ) = 

for k = 1,...,N— 1. The resulting update scheme is not self-starting. To initialize the 
method one needs to provide (</o,9i) € Q x Q as opposed to a point (qo,Vo) € TQ. 

Symplecticity As in the continuous theory symplecticity follows from restricting <& d (oj, •) 
to pathwise unique solutions of the stochastic discrete Euler-Lagrange equations, <& d . Since 
pathwisc unique solutions can be parametrized by initial conditions we regard the restricted 
action as <3 d : ft x Q x Q — ► R. Taking its first variation gives, 

N-l 

d<S d (u,q ,qi) ' {SQo,6Qi) = X! ( D 3 L d(tk,u,q k ,q k+ i) + D 4 L d (t k - 1 ,uj,q k - 1 ,q k )) ■ Sq k 

fe=i 

+ D 3 L d (t ,uj,qo,qi) ■ 5q + D 4 L d (t N -\,u,q N -i,q N ) ■ Sq N . 

Because of the restriction to solution space the sum vanishes, and the boundary terms 
remain: 

d<S d (u), q , gi) • (Sq , Sq 4 ) = D 3 L d (t , w, q , qi) ■ Sq + D 4 L d (t N -i,uj, qN-i, 9jv) ■ Sq N . 

These boundary terms define left and right one-forms as follows: 

<d + (t k ,uj,q k ,q k+1 ) ■ (Sq k ,5q k+1 ) = D 4 L d (t k ,uj, q k , q k+1 ) ■ Sq k+1 
@~(t k ,u),q k ,qk+i) ■ (5q k ,6q k+1 ) = D 3 L d (t k ,u),q k ,q k +i) ■ Sq k 
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which from d 2 Ld = satisfy: 

de+ = -de - = n. 

Applying the second exterior derivative to 25 ^ implies, 

fi(tjV-l, CJ, qN-l,qN){^N-l' S ^N)i S <lN-l^ S QN) = ^(<o> ^, 9o, (<5<7o . S 1l ) (<^0 . ) 

since d 2 ©^ = 0. Hence, the discrete flow preserves the symplectic form fi. 

Discrete Momentum Map Consider the left action of a Lie Group G on Q. If the 
stochastic discrete Lagrangian is infinitesimally symmetric, then the associated momentum 
map is preserved. A sufficient condition for this is that the discrete Lagrangian is invariant 
with respect to the left action of G. The proof is sketched out here. 

Let the action on the discrete configuration manifold be denoted by &q x q ■ GxQxQ — > 
Q x Q and defined by: 

^QxQ(g,qi,q2) = (®(g,qi),$(g,q2))- 
The associated infinitesimal generator is denoted £q x q : Q X Q —* T(Q x Q) and by definition 

£.QxQ{<h.,(&) = j- *QxQ(exp(s^),gi,g 2 )| s=0 ■ 
Assume that Ld is infinitesimally symmetric, i.e., 

dL d ■ £,q x q = 6+ • Cqxq + 0~ • Cqxq = 0. 
By this condition the left and right discrete momentum maps, J + , J~ : Q x Q — > g*,: 

j- ■ z = -e- • Cqxq 

are equal, i.e., J + = J - = J. Consider the restricted action sum and compute its differential 
in the direction of the infinitesimal generator to obtain: 

d<5 d (w, qo,qi) ■ ^QxQ(qo,qi) = e~(£ ,^,go,<7i) • ^QxQ(qo,qi) 

+ Q + (t N _ 1 ,LU, qjv-l, <ZJv) ' £,QxQ(qN-l,qN), 

which can be rewritten in terms of the momentum maps as 

d&d(uj,qo,qi) ■ £qxq(<zo,<?i) = -J~{to,w,qo,qi) ■ £, + J + {t N -i,^,qN^i,qN) ■ £, = 0. 

Since the left and right momentum maps evaluated at the same point are equal, the mo- 
mentum map J is preserved under the discrete flow. 

5 Stochastic Variational Euler Integrator 

In the deterministic setting, the HP context provides a practical way to design discrete 
Lagrangians and obtain one-step methods on TQ or T*Q as pointed out in [3]. In this 
section we examine variational Euler methods extended to the stochastic context following 
the continuous stochastic HP theory laid out in §2. 
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Stochastic Variational Euler on W 1 To discretize the stochastic HP action integral one 
needs to replace the continuous Lagrangian, stochastic integral and kinematic contraint by 
discrete approximants. We begin by introducing a first-order discretization of the kinematic 
constraint in (3.1). Let [a,b] and N be given and define the fixed step size h = (b — a)/N 
and tk = hk, k = 0, N. 

A discretization of the kinematic constraint can be obtained by introducing a discrete 
sequence {q k } k=0 sucn that € Q and a finite difference map ip : Q x Q — > TQ. For 
example, if Q is a vector space the following backward difference map can be introduced: 



<p(qk,qk+i) = ( qk 



qk+i - qk 



1, ui and fc 



Let ~ A/"(0, /i) be normally distributed random variables for ; 
0, ...,N — 1. In terms of the discretization of the kinematic constraint, the corresponding 
discrete HP action sum takes the following simple form: 



N-l 



«5=E 



C(q k ,v k )h + y^7i(gfc)-Bf + (Pfc+i, (?fc+i - - Ufc+i) h 



k=0 L t=l 

The stochastic discrete HP equations are given by: 
q k+ i = q k + hvk+i, 



Pfc+i =p k + h — (q k ,v k ) + 2^ -Q-{qk)Bi , 



<9g 

Pit = -^{qk^kj- 



Stochastic Variational Euler on Lie Groups In the context of Lie groups, the recon- 
struction equation is discretized using canonical coordinates of the first kind, r : q — ► G, 
as explained in [3; 4]. As in the vector space case, we define a finite difference map 
(p : G x G — > G x g that provides a first-order approximant to the reconstruction equa- 
tion: 

<f{9k,gk+i) = (gk+i,T- 1 (g~ k 1 g k +i)/h) e G x g. 
A first order approximant to the stochastic left-trivialized action integral is given by, 



N-l 

E 

k=0 



K9k-,(,k)h + ^"i l (g k )B i [ + (nk+i,r 1 (g k 1 g k +i)/h - £ fc+1 ) h 



i=l 



(5.1) 



Let dr 1 : g x g — > g denote the right trivialized tangent of t 1 as defined in [4]. The 
stochastic left-trivialized discrete HP equations are (cf. definition 3.4): 

gk+i = g k T(h£ k+1 ), 

m 

( dT hl +1 T^+i = (d^JVfc + hl g {g k ,i k ) + 5^(7i) fl (Sfc)B*, 

f=i 
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Holonomic Constraints Si Nonconservative Effects Holonomic constraints can be 
added via discrete Lagrange multipliers, and nonconservative effects via discrete impulses 
as described below in the Lie group context. Suppose that G is an n-manifold and that 
the mechanical system evolves on a submanifold S C G defined as the zero level-set of 
(p : G — > M. k where k < n and S — y> _1 (0). Further suppose that there exists a force 
field F : G — > T*G. These effects are appended by using the following modified stochastic 
left-trivialized HP action principle: 

JV-l r m 

5 l ^k, Zk)h + li{9k)B^ + (nk+uT-^g^gk+J/h - h + (X k ^(g k )) h 

k=0 L i=l 
N-l 

+ F (9k) ■ Sg k h = 0. 

fc=0 

This algorithm for G = M. n is the stochastic analog of constrained symplectic Euler, and the 
numerical analysis of this method is discussed in [6]. 



6 Langevin-type equations for Multiple Rigid Bodies 

Continuous Description Consider a mechanical system consisting of K rigid bodies 
interacting via a pairwise potential dependent on their positions and orientations. Let 
(xj(f), Vi(t), Ri(t), u>i(t)) S TSE(3) denote the translational position, translational velocity, 
orientation, and spatial angular velocity of body i where i = 1, K. Let rrii and I; denote 
the mass and diagonal inertia tensor of body i. The left-trivialized Lagrangian for the 
system is given by: 



K ! 

l(xi,Vi,Ri,uji) = ^ ~Y v 'i v i + RihRj i^t - U(xi,Ri). 



i=l 

Note that l(xi, Vi,Ri,u>i) is shorthand notation for l(x\, V\, R\,ui\, ■ ■ ■ , Vk, Rk,0Jk)- We 
will use this shorthand notation elsewhere to simplify the expressions. The path that the 
stochastic mechanical system takes on the time-interval [a, b] is one that extremizes the HP 
action: 

m / dx \ I dR \ 

l{xi,Vi,Ri,u>i)dt + y]7q(xj,Ri) o dW q + yjj, - Vij dt + "^-^j 1 " dt 



for arbitrary variations with fixed endpoints: (xi(a), Ri(a)) and (xj(6), Ri(b)). The corre- 
sponding SDEs of motion are given by: 

dxi = Vidt (reconstruction equation), 

tn 

dpi = —U Xi dt + ^(7^)2;; dW q (t, u>) (stochastic EL equations), 

q=l 

Pi = rriiVi (Legendre transform), 

dRi = QiRidt (reconstruction equation), 

m 

dir.i = —U^dt + V^(7g)ij 4 dWq (stochastic LP equations), 

q=l 

Hi = RiliRjuJi (reduced Legendre transform). 
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for i=l, K . The terms U Xi and [/r ; arc defined in terms of the inner product on R 3 as, 

Uj t y = (j^.iV) = d Xi U(xi,Ri) ■ y, 

where d R V : S0(3) -> T£. S0(3), and d^t/ : M 3 -► T^M 3 as defined below. Adding 
dissipation so that the Gibbs distribution is invariant under the stochastic process defined 
by the above SDE with dissipative drift, yields Langevin-type equations for rigid-body 
systems (see, e.g., [5]). 

Stochastic Variational Euler integrator For the discrete description, the variational 
Euler integrator provided earlier is implemented. Let ~ Af(0, h) be normally distributed 
random variables for q — 1, m and k = 0, N — 1. The action sum is given by 

N-l 



fc=0 

Tfc+i __— 1 z' r)fe-(-l / r?fe\T\ /l , ,fc+l 



9=1 



Stationarity of this action sum implies the following discrete scheme: 

= a * + hv k ; +1 

m 

P- +1 =v\- hU Xt (a:*, R>°) + R*)B k q 

9=1 

Pi = w 



= r R* 

* * m 



fc 

9 

9=1 

fc _ ofcir / ofc^T, ,fc 



for i = 1, ...,K. Assuming the Legendre transforms are invertible, this integrator has the 
attractive property that the translational and rotational configuration updates, and the 
translational momentum update, are explicit. One only has to perform an implicit solve for 
the discrete Lie-Poisson part. Even that computation is straightforward since the torque due 
to the potential is only a function of the orientation and position at the previous time-step. 
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